simulation results/choosing_k.R

In BriceonWiley/IntegratedLikelihood.R:

map(
  seq(50, 700, 25),
  ~ simulate_p(
    N = .x, p = 0.1, phi = 0.1, theta = 0.1, k = seq(0.05, 0.65, 0.05)
  )
) -> new_sim_res_k

mix <- expand_grid(N = seq(50, 700, 50), p = seq(0.1, 0.9, 0.1))
mcMap(
  function(N, p) {
    simulate_p(
      N = N, p = p, phi = 0.1, theta = 0.1, k = seq(0.05, 0.65, 0.05)
    )
  }, N = mix$N, p = mix$p, mc.cores = 3
) -> table_results
mcMap(
  function(N, p) {
    simulate_p(
      N = N, p = p, phi = 0.5, theta = 0.5, k = seq(0.05, 0.65, 0.05)
    )
  }, N = mix$N, p = mix$p, mc.cores = 3
) -> table_results_050

coverage_summary <- map_dfr(table_results_050, ~ coverage_p(.x, TRUE))

saveRDS(table_results_050, "/Users/briceon_wiley/Documents/RStudio/comparingCIs.R/simulation results/results_p_adj/big_results_050.rds")
saveRDS(coverage_summary, "/Users/briceon_wiley/Documents/RStudio/comparingCIs.R/simulation results/results_p_adj/summarized_big_results_050.rds")

choose_best <- function(data) {
  data %>%
    mutate(distance = abs(Coverage - 0.95)) %>%
    filter(distance == min(distance)) %>%
    select(-distance)
}

coverage_summary %>%
  select(-mW) %>%
  gather("Type", "Coverage", -N, -p) %>%
  group_split(N, p) %>%
  map_dfr(choose_best) -> new_results

list(
  "choice_of_k" = new_results %>%
    select(-Coverage) %>%
    mutate(
      Type = case_when(
        (Type == "ILR") ~ "1.00",
        TRUE ~ str_sub(Type, 3L, 6L)
      )
    ) %>%
    spread(p, Type),
  "coverage" = new_results %>%
    select(-Type) %>%
    spread(p, Coverage)
) -> k_summary

saveRDS(k_summary, "/Users/briceon_wiley/Documents/RStudio/comparingCIs.R/simulation results/results_p_adj/k_050.rds")